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Abstract - Bone is a biomaterial undergoing continuous renewal. The renewal process is known as bone 
remodelling and is operated by bone-resorbing cells (osteoclasts) and bone-forming cells (osteoblasts). An 
important function of bone remodelling is the repair of microcracks accumulating in the bone matrix due 
to mechanical loading. Cell-cell communication between cells of the osteoclastic lineage and cells of the 
osteoblastic lineage is thought to couple resorption and formation so as to preserve bone integrity and achieve 
homeostatic bone renewal. Both biochemical and biomechanical regulatory mechanisms have been identified in 
this coupling. Many bone pathologies are associated with an alteration of bone cell interactions and a consequent 
disruption of bone homeostasis. In osteoporosis, for example, this disruption leads to long-term bone loss and 
fragility, and can ultimately result in fractures. 

Here we focus on an additional and poorly understood potential regulatory mechanism of bone cells, that 
involves the morphology of the microstructure of bone. Bone cells can only remove and replace bone at a 
bone surface. However, the microscopic availability of bone surface depends in turn on the ever-changing 
bone microstructure. The importance of this geometrical dependence is unknown and difficult to quantify 
experimentally. Therefore, we develop a sophisticated mathematical model of bone cell interactions that 
takes into account biochemical, biomechanical and geometrical regulations. We then investigate numerically 
the influence of bone surface availability in bone remodelling within a representative bone tissue sample. 
Biochemical regulations included in the model involve signalling molecules of the receptor-activator nuclear 
factor kB pathway (rank-rankl-opg), macrophage colony-stimulating factor (mcsf), transforming growth 
factor p (tgf|3) and parathyroid hormone (pth) . For the biomechanical regulation of bone cells, a multiscale 
homogenisation scheme is used to determine the microscopic strains generated at the level of the extravascular 
matrix hosting the osteocytes by macroscopic loading. The interdependence between the bone cells' activity, 
which modifies the bone microstructure, and changes in the microscopic bone surface availability, which in turn 
influences bone cell development and activity, is implemented using a remarkable experimental relationship 
between bone specific surface and bone porosity. Our model suggests that geometrical regulation of the 
activation of new remodelling events could have a significant effect on bone porosity and bone stiffness. On 
the other hand, geometrical regulation of late stages of osteoblast and osteoclast differentiation seems less 
significant. We conclude that the development of osteoporosis is probably accelerated by this geometrical 
regulation in cortical bone, but probably slowed down in trabecular bone. 

Key words: mechanical feedback, geometrical feedback, specific surface, bone remodelling, bone stiffness, 
osteoporosis 



1 Introduction 

Bone is a biomaterial that has a variety of physiological 
functions. In addition to load bearing and support for 
locomotion, bone protects internal organs and participates 
in calcium and phosphorous homeostasis. From an en- 
gineering perspective the structural function of bone is 
most importantly characterised by its stiffness and strength. 
Daily activities (such as walking and running) subject bone 
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to periodical loads which, over extended periods of time 
(weeks, months and years), can lead to fatigue damage and 
the formation of microcracks. If these microcracks are not 
removed in due time, their evolution may result in a macro- 
scopic structural failure, i.e., a fragility fracture. To prevent 
the occurrence of fatigue fractures, nature has equipped 
bone tissues with a cellular mechanism of self-repair [1], 
referred to by biologists as 'bone remodelling' [2,3]. Bone 
remodelling is a coordinated process of bone resorption by 
cells called osteoclasts, and bone formation by cells called 
osteoblasts. Osteoclasts and osteoblasts usually operate 
together in self-contained groups processing the renewal 
of a localised portion of the bone tissue. These groups 
are called bone multicellular units (bmus) and constitute 



a single 'remodelling event'. There are about 1.7-10 such 
bmus in a normal adult skeleton [2-4]. Cell population 
and cell activity in a bmu are tightly controlled to establish 
local bone homeostasis (i.e., balanced bone resorption and 
bone formation) . In bone pathologies, this cellular control 
is perturbed and homeostatic bone renewal is disrupted. 
In osteoporosis, bone is progressively lost, which results in 
reduced bone stiffness and strength. 

Over the last decades, bone biologists have identified 
a large number of biochemical regulatory factors influ- 
encing bone remodelling. The formation of osteoclasts 
has been shown to rely crucially on macrophage colony- 
stimulating factor (mcsf) and on the receptor- activator 
nuclear factor kB (rank) cell signalling pathway, which 
involves the receptor rank, the ligand rankl and osteo- 
protegerin (opg) [5,6]. Rankl activates the rank receptor 
on precursor osteoclasts, which triggers their development 
and sustains their activity The soluble molecule opg is 
a decoy receptor of rankl which can prevent rankl from 
binding to rank. Another important molecule mediating 
the communication between osteoblasts and osteoclasts is 
transforming growth factor |3 (tgf|3). Tgfp is stored in 
high concentrations in the bone matrix. It is released into 
the bone microenvironment, where it exerts its action on 
several bone cells, during bone matrix resorption by active 
osteoclasts [6] . The existence of a mechanical regulation of 
bone remodelling has long been suspected. It is now well 
established that mechanical feedback is a key regulatory 
mechanism to maintain bone mass [7-12]. The commonly 
accepted view is that osteocytes act as mechanosensors 
that transduce local mechanical signals into biochemical 
responses. These biochemical responses are thought to 
regulate the initiation of bone remodelling processes and 
to modulate the coupling between bone resorption and 
formation (see e.g. [1] and Refs cited therein). 

The existence of biochemical and biomechanical reg- 
ulations of bone cells is well-established and has been 
extensively studied. However, the notion that the mor- 
phology of the microstructure of bone may induce an 
additional regulation of bone cells of purely geometrical 
nature is not often mentioned in the recent literature. This 
may be due to the experimental difficulty of assessing the 
importance of a geometrical regulation. Biochemical and 
biomechanical regulations can experimentally be partially 
or fully repressed by selective gene knock-outs or mono- 
clonal antibodies targeting key components in the bone 
cell signalling pathways. By contrast, one cannot simply 
"switch off" a geometrical regulation of bone remodelling 
when this self-repair process modifies the microstructure 
(and so the geometry) of the material. 

Bone tissue is diverse and exhibits a broad variety of 
microstructures. However, two distinctive types of bone 
tissue are usually identified: cortical bone and trabecular 
bone [3] (see images in Figure 1). Cortical bone has typical 
porosities of 0.05-0.15 while trabecular bone has typical 
porosities of 0.65-0.85 [2,3]. Mathematical models for 
the estimation of mechanical properties of bone tissue have 
shown that bone stiffness is predominantly determined by 
the porosity f^, 1 the interaction of the different mate- 

1 The total porosity of bone is made of a vascular porosity, which con- 



rial phases and pore shape, while other microstructural 
characteristics such as the exact pore distribution play a 
secondary role [14-16]. For biochemical processes, pore 
morphology can be expected to play a significant role. 
Indeed, pore morphology determines the so-called specific 
surface S v (i.e., the amount of bone surface available in 
a representative volume element), which is an essential 
geometrical factor for the bone cells. Bone cells require 
a bone surface to fulfill their functions, whether to initiate 
a bone remodelling process or to operate resorption and 
formation. Osteoclasts require attachment to a particular 
area of the bone surface before resorbing. Osteoblasts are 
observed to only secrete osteoid (a collagen-rich substance 
which later mineralises and becomes new bone matrix) on 
existing bone surfaces. Finally, mechanical signals sensed 
by osteocytes embedded in the bone matrix are passed 
on to bone cells in the vascular cavity through the bone 
surface. Effects similar to chemical exchange reactions 
between pore walls and solutes in fluid-saturated porous 
materials can be expected to occur in this context. 

The issue of quantifiying the role of bone surface avail- 
ability in bone remodelling was raised by bone biologists 
already some time ago [2, 13, 17]. In Ref. [13], Martin 
provides a first attempt to investigate theoretically the 
effect of a geometrical regulation of bone remodelling in 
osteoporosis (see Figure 1). Osteoporosis is associated 
with increased porosity in both cortical and trabecular 
bone [18]. In Martin's own words: "In [cortical] bone, 
increased porosity provides more surface area on which cells 
can work, thereby increasing the capacity for further porosity 
changes. In [trabecular] bone, increased porosity decreases 
the amount of surface available to the cells, thereby decreas- 
ing the capacity for further remodelling." 

While the proposed mechanism of geometrical feedback 
on bone remodelling seems plausible, it is difficult to test 
its validity experimentally and to determine its impor- 
tance quantitatively. Some researchers have employed the 
concept of geometrical feedback for simulations of bone 
remodelling [21]. However, to our knowledge, there is no 
systematic study in the literature of the effects of a possible 
geometrical regulation at several stages of the remodelling 
sequence. Also, the interplay between geometrical feed- 
back and mechanical feedback in bone remodelling has not 
been investigated. A mechanical feedback has the potential 
to stabilise bone loss or gain [7, 22] and may either com- 
pete with or enhance the effect of the geometrical feedback 
seen in Figure 1 depending on the type of bone. 

The aim of this paper is to address the above questions 
using a state-of-the-art computational model of bone re- 
modelling. The following questions related to geometrical 
feedback are investigated: (i) At which stage of the bone 
remodelling sequence (i.e., activation, resorption, forma- 
tion) does geometrical feedback have the strongest effect?; 

tains marrow components, blood vessels, bone cells and their precurors, 
and the lacunae-canaliculi porosity, which contains osteocytes and their 
processes. The lacunae-canaliculi porosity is only a small fraction of the 
total porosity (see e.g. [13, Table 1]) and no remodelling occurs at these 
surfaces. Therefore, the lacunae-canaliculi porosity will not be consider 
in this work and we will refer to the vascular porosity simply as the bone 
porosity. Similarly, in the present context, we are not interested in the 
intercrystalline and intermolecular porosities, which we simply regard as 
part of the 'solid bone matrix'. 
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Figure 1 - Schematic representation of a possible effect of geometrical 
feedback on the evolution of vascular porosity (/ vas ) in osteoporosis 
both in cortical bone (lower curves) and trabecular bone (upper curves) 
according to Martin [13]. Dashed curves show the linear increase 
in porosity that is obtained without geometrical feedback, while the 
solid curves incorporate geometrical regulation (a constant pathological 
skeletal imbalance of -2 p.m/yr is assumed). Typical bone microstructures 
for cortical bone (bottom) and trabecular bone (top) in normal subjects 
(left) and osteoporotic subjects (right) are also represented. 

(ii) How do geometrical and mechanical feedbacks inter- 
act?; (iii) What is the impact of geometrical feedback in 
osteoporosis in terms of bone porosity and bone stiffness? 

To address these questions we extend a previously de- 
veloped mathematical model of bone remodelling [22-24]. 
This multiscale model takes into account both biochemical 
and biomechanical regulations of bone remodelling. Bio- 
chemical regulatory factors include the rank-rankl-opg 
pathways together with the action of tgf|3 on bone cells. 
Biomechanical regulation of bone formation and bone 
resorption is mediated by the microscopic strain energy 
density (sed) of the bone matrix. This strain energy den- 
sity is calculated from a micromechanical homogenisation 
scheme. The introduction of geometrical feedback due to 
microscopic bone surface availability is elucidated through 
a phenomenological relationship between the specific sur- 
face and the vascular porosity obtained from various types 
of bone [13]. 

2 Mathematical model of bone re- 
modelling 

Few mathematical models of bone remodelling include 
explicitly biochemical interactions of bone cells that couple 
bone resorption and bone formation. Lemaire etal. [25] 
have proposed a bone cell population model that incorpo- 
rates some of the most important known bone biology. This 
cell population model was further developed by Pivonka 
etal. to investigate the effect of rankl expression on os- 
teoblasts of varying maturity [23]. These models have 
been shown to correctly capture important physiological 
behaviours of bone remodelling both in bone homeostasis 
and in bone pathologies [24, 26-28] . Recently, we have 
proposed an extension of the model of Ref. [23] to include 
a biomechanical regulatory mechanism to trigger bone cell 
responses, leading to a fully coupled model of biochemical 
and biomechanical regulations [22]. This model uses 



a novel multiscale approach based on micromechanical 
homogenisation of bone stiffness. It allows to consistently 
calculate bone matrix strains at the osteocyte level. Osteo- 
cytes convert this micromechanical signal into biochemical 
signals to bone cells in the vascular space, effectively 
providing a biomechanical feedback regulation of bone 
remodelling. 

In the following, we present an extension of the model 
of Ref. [22] that incorporates the influence of microscopic 
bone surface availability at various stages of the bone 
remodelling sequence. 

2.1 Fundamental biochemical and biome- 
chanical regulation of bone remodelling 

The biochemical regulatory mechanisms considered in the 
model have been described in detail in Refs. [23, 24]. 
Three stages of osteoblast development are considered: 
uncommitted osteoblast precursors (ob u s), pre-osteoblasts 
(oB p s) and active osteoblasts (oB a s). Similarly, three stages 
of osteoclast development are considered: uncommitted 
osteoclast precursors (oc u s), pre-osteoclasts (oc p s) and 
active osteoclasts (oc a s). Figure 2 schematically depicts 
these bone cell types along with their biochemical, biome- 
chanical and geometrical regulations. These regulations 
are summarised sequentially in the following. 

The generation of osteoblasts is assumed to be regulated 
by transforming growth factor |3 (tgf|3) (a factor released 
from the bone matrix during osteoclastic bone resorption) 
while the generation of osteoclasts is assumed to be regu- 
lated by rankl and opg (molecules in the receptor-activator 
nuclear factor kB (rank) system, that are expressed by 
osteoblasts). Tgf(3 promotes the differentiation of uncom- 
mitted osteoblast progenitors (ob u ) into pre-osteoblasts 
(oB p ), but it inhibits the differentiation of pre-osteoblasts 
into active osteoblasts (0B a ). Furthermore, tgf|3 promotes 
osteoclast apoptosis (programmed cell death). Rankl is a 
protein expressed on the surface of pre-osteoblasts. The 
binding of rankl to the receptor rank found on pre- 
osteoclasts promotes the differentiation of pre-osteoclasts 
into active osteoclasts. However, this binding may be 
prevented by the presence of opg, a decoy receptor of 
rankl produced in soluble form by active osteoblasts. Fur- 
thermore, the circulating parathyroid hormone pth induces 
rankl expression and downregulates opg expression by 
osteoblasts to produce a systemic increase in available 
rankl, resulting in an increase in the formation and ac- 
tivity of osteoclasts. The differentiation of uncommitted 
osteoclast progenitors (oc u s) into pre-osteoclasts (oc p s) 
requires signalling by both rankl and macrophage colony- 
stimulating factor (mcsf). All these signalling actions are 
accounted for using mass action kinetics for the chemical 
bindings between ligands and their receptors [29]. The 
fraction of occupied receptors on a cell (i.e., bound to a 
ligand) is assumed to determine the strength of the signal 
received by the cell, and so is assumed to determine the 
strength of the cell's response. 

The biomechanical regulatory mechanisms considered in 
the model are described in detail in Ref. [22]. Mechanical 
disuse is known to increase bone resorption by increasing 
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Figure 2 - Cell population model of bone remodelling, including several developmental stages of osteblasts (ob u , 



, OB a ) and osteoclasts (oc u , 



oc a ) and their biochemical regulation (pth, rank-rankl-opg, tgf|3, mcsf), biomechanical regulation (sed, £b m ) and geometrical regulation (S v ) (see 
text for further explanations) . 



the rankl/opg ratio, which increases the differentiation 
of pre-osteoclasts into active osteoclasts. In the model, 
this is implemented using a mechanically controlled rankl- 
production term, see Eq. (10). Mechanical overuse is 
believed to increase bone formation by stimulating wnt 
signalling to pre-osteoblasts, which increases their prolif- 
eration and ultimately leads to an increased population 
of active osteoblasts [22,30]. In the model, this is im- 
plemented using a mechanically controlled proliferation of 
pre-osteoblasts, see Eq. (11). 

2.2 Geometrical and morphological charac- 
teristics of bone: porosity and specific 
surface 

To represent the microstructure of bone at the tissue level, 
the most important geometrical and morphological param- 
eters are the vascular porosity (f vas ) and the specific surface 
(S v ). In cortical bone, vascular porosity corresponds to 
the so-called 'Haversian canal' system. In trabecular bone, 
vascular porosity corresponds to the marrow space around 
the trabecular struts [3]. The vascular compartment con- 
tains all the bone cells considered in our model. Vascular 
porosity is defined as the volume fraction of vascular pores, 
i.e., the volume of vascular pores (V* vas ) per tissue volume 
(V r ): 2 



/vas — Kas/^T- 



(1) 



The bone matrix volume fraction is defined in the same way 
as the volume of bone matrix (Y hm ) P er tissue volume, i.e.: 



fbm — V bm/ V T- 



(2) 



We recall that all porosities at observation scales below 
the vascular porosity, such as the lacunar porosity, and 
the canaliculi connecting the lacunae, are not involved 
in the remodelling process, i.e. the intricate biochemi- 
cal processes take place within the vascular pore space 



2 The tissue volume V T is assumed to be of the order of 1-3 mm 3 . 
This volume corresponds to a representative volume element (RVE) large 
enough to comprise several remodelling events (bmus) [2], yet small 
enough to represent spatial heterogeneity of bone tissue (in particular, 
small enough to distinguish cortical bone and trabecular bone) [16]. 



(see also footnote 1). From Eqs (l)-(2), it follows that 
/vas +/bm = 1- Typically cortical bone exhibits a range of 
porosities / vas 0.05-0.15 while trabecular bone exhibits 
a range of porosities / vas 0.65-0.85. 

The specific surface (or surface density) of a porous 
material is defined as the interstitial surface area of the 
pores (S p ) per tissue volume, i.e.: 



SJV T 



(3) 



with dimensions [mm 2 /mm 3 ]. Generally speaking, the 
specific surface is an important quantity for a variety of 
phenomena in porous media. For example, it determines 
the adsorption capacity of industrial adsorbents and plays 
an important role in determining the effectiveness of cat- 
alysts and ion exchange columns and filters. It is also 
related to the fluid conductivity or permeability of porous 
media (see e.g. Ref. [31]). Experimentally, S v is commonly 
estimated by adsorption methods, quantitative stereology 3 
fluid flow and micro-computed tomography. As mentioned 
above, in bone the specific surface is important as it 
determines the available working area for osteoblasts and 
osteoclats. The specific surface can also be expected to have 
an influence on the transmission of specific signalling by 
osteocytes in the bone matrix to osteoblasts and osteoclasts 
developing in the vascular pores. 

The microstructure of a material determines both the 
porosity and the specific surface. Depending on the mi- 
crostructure, different materials exhibit different values 
for these quantities. Bone tissue spans a wide range of 
porosities, each of which is characteristic of a particular 
micro-architecture, and so of a particular value of the 
specific surface. Based on a large number of experimental 
data, Martin has provided a remarkable phenomenological 
relationship between bone specific surface (S v ) and vascu- 
lar porosity (/ vas ) [13, Eq. (68)]: 



SvC/vas) — a /vas + ^ fvas" + c /vas 3 + d f vas 4 + e /v: 



(4) 



where the polynomial coefficients are estimated as a = 
32.3, b = -93.9, c = 134, d = -101, and e = 28.8 (in 
[mm -1 ]). In Figure 3, the relation (4) is plotted together 



3 Stereology is the study of three-dimensional properties of objects 
observed in two-dimensional sections. 
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Figure 3 - Relation between bone specific surface (_S V ) and vascular 
porosity (/ vas ) (modified from Martin [13, Fig. 19]). Data points were 
measured on histological sections and microradiographs of different types 
of bone. Solid symbols: healthy bone; open symbols: diseased bone 
(osteoporosis, osteogenesis imperfecta, and osteomalacia); circles: human 
femoral bone; squares: human iliac crest; diamonds: human rib; trian- 
gles: human vertebra. The maximum of the curve S v (f vas ) (thick solid 
line, Eq. (4)) at / vas = / v * ss 0.37 is indicated by dotted lines. The dashed 
curves are obtained by varying the polynomial coefficients in Eq. (4) to 
enclose most of the range of measurements. 

with experimental data obtained from various types of 
human bone (femur, iliac crest, vertebra, rib) both in health 
and disease. The data and the curve S v (/ vas ) show two 
important characteristics: (i) all specimens approximately 
follow the same S v (/ vas ) curve independently of bone type 
and with no significant difference between diseased and 
healthy bone. This remarkable universality establishes the 
curve S v [f vas ) as an intrinsic property of bone; (ii) the 
specific surface exhibits a maximum at a porosity /* of 
about 0.37, intermediate between cortical and trabecular 
bone (bone at such porosity is denoted as 'porous-compact' 
bone). 

In this work, we use Eq. (4) to include geometrical 
regulation in the model as follows. The evolution of the 
bone cell populations predicts the evolution of the vascular 
porosity due to resorption by osteoclasts and formation by 
osteoblasts. Eq. (4) then enables us to estimate changes in 
the specific surface associated to the changes in porosity. 
This change in microscopic bone surface availability is 
in turn assumed to influence the evolution of the cell 
populations. 

2.3 Bone cell governing equations 

In the model, osteoblasts and osteoclasts of different de- 
velopmental stages are considered. The populations of os- 
teoblasts and osteoclasts are therefore heterogeneous. The 
composition of these populations is determined by follow- 
ing individually the populations of each of the developmen- 
tal stages of osteoblasts and osteoclasts mentioned above. 
The evolution of these cell (sub) populations is transcribed 
mathematically as so-called 'rate equations' [22,23]. The 
populations of uncommitted progenitor osteoblasts (ob u s) 
and uncommitted progenitor osteoclasts (oc u s) are as- 
sumed constant and so are not state variables. These 
uncommitted cells represent a pool of progenitor (or stem) 
cells that is assumed to be maintained by self-renewal 



unlimitedly. The possibility for geometrical regulation is 
included at each stage of osteoblast and osteoclast devel- 
opment through functions of the specific surface, namely 
g B u CS v ), g B v (.Sv)> goc u (S v ) and &oc p CS v ) (see Figure 2). In 
the following, we denote the bone cell densities within a 
tissue sample (number of cells per unit volume) by their 
symbol ob u , ob p , OB a , oc u , oc p , oc a . The concentrations of 
the biochemical signalling molecules within a tissue sample 
(number of molecules per unit volume) is also denoted 
by their symbol tgf(5, rankl, etc. 4 Based on the above 
descriptions of the biochemical, biomechanical, and geo- 
metrical regulatory mechanisms, the governing equations 
of the bone cell densities in the model are expressed as: 

^OB p = {g OBu (S v )D OBu 7l^ 0Bu } OB u + {P OBp n a t% p } OB p 

-{go Bp (S v )D OBp < G ;f OBp }oB p (5) 
^OB a = {g 0Bp (S v )D 0Bp 7i™ F f 0Bp } OBp -A QBa OB a (6) 

^oc p = {g 0Cu (S,)D 0Cu Cu n™^} oc u 

-{goc p (Sv)0 OCp <™}oc p (7) 

A oc = S g (S V )D 0C 7T RA ™ L } OC D - \A 0C 7I Tt f } OC a . 
^£ a ^e>oc p \. vj oCp act,oc p J P ^ oc a act,oc a J a 

(8) 

Below, we discuss the various quantities occurring in these 
equations (see Refs [22, 23] for more details). The pa- 
rameters D OB , D 0B , D oc , D oc denote differentiation rate 
parameters for uncommitted osteoblast progenitors, pre- 
osteoblasts, uncommitted osteoclast progenitors and pre- 
osteoclasts, respectively; P OB is a proliferation rate param- 
eter for pre-osteoblasts; A 0Ba and A oc ^ are apoptosis rate 
parameters for active osteoblasts and active osteoclasts. 
Biochemical regulation of cell development is achieved 
through so-called "activator" and "repressor" functions 
of the biochemical signalling molecules tgfp, mcsf and 
rankl as follows. The functions 7i T ^ 0B ^, 71 rep™ ' an d 

^act.oc are activator and repressor functions regulating 
osteoblast differentiation and osteoclast apoptosis based 
on the concentration of tgfr. The function n MC f S¥ is an 

1 act,oc u 

activator function regulating differentiation of oc u s into 
oc p s based on the concentration of macrophage colony 
stimulation factor (mcsf). The functions 7i™^ and ti™!; 

aCt,OC L | 3.Ct,OCp 

are activator functions regulating osteoclast differentiation 
based on the concentration of rankl. For simplicity, we 
assume that 7r™ = 7rjJ£^ • The equations governing 
the evolution of the biochemical signalling molecules tgf|3, 
mcsf, rankl, rank, opg and pth and the form of the activa- 
tor and repressor functions are presented in Appendix A, 
see Eqs (15)-(22). Tables 1-3 in Appendix A list the 
dynamic quantities, the biochemical parameters and the 
biomechanical quantities of the model. 



4 To align with common practice, we shall use the terminology 'density' 
for cells and 'concentration' for signalling molecules, even if the units are 
chosen the same. 
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Figure 4 - Microscopic strain energy density tp^ m derived from the 
micromechanical homogenisation procedure as a function of vascular 
porosity / vas for the case of a uniaxial compressive loading £ = 2 um e 3 ® e 3 
with S uni = -30 MPa. 



2.3.1 Mechanical feedback regulation 

Biomechanical studies suggest that the strain energy den- 
sity is an important quantity that determines bone adap- 
tation to various mechanical loads. This quantity is com- 
monly chosen in the literature due its scalar nature as a 
measure of mechanical stimulus sensed by the bone cells to 
drive bone adaptation [3,32]. 

The model of biomechanical regulation developed in 
Ref. [22] similarly uses the strain energy density as the 
signal conveying mechanical information to the bone cells. 
However, as discussed in Ref. [22], it is paramount to 
estimate consistently the bone matrix strain energy density 
at the micro-scale where osteocytes sense this mechanical 
signal and transduce it into a biochemical response. In 
Ref. [22], we used a homogenisation procedure based 
on Eshelby's classical matrix inclusion problem and the 
Mori-Tanaka scheme to estimate the microscopic strains 
e bm generated at the level of the extravascular matrix 
hosting the osteocytes, by the macroscopic loading of the 
tissue. This homogenisation procedure leads to an intricate 
dependence of the microscopic strain energy density of the 
bone matrix, ipbm> upon the macroscopic stress tensor S 
and the vascular porosity / vas : 



V'bm = ^bm(S,/ V as)- 



(9) 



In this paper, we will exemplify the effect of biomechanical 
regulation for the situation of a constant uniaxial com- 
pressive loading only (i.e., S = S um e 3 ® e 3 with E um = 
—30 MPa). For this situation, the dependence of ip hm upon 
/ vas is plotted in Figure 4. The reader is referred to Refs [14, 
15, 22, 33] for the general mathematical specification of 
i^ bm and a detailed presentation of the homogenisation 
procedure. As mentioned in Section 2.1, the biomechan- 
ical regulation of bone remodelling is believed to operate 
through different pathways for resorption and formation. 
Following Ref. [22], the biomechanical regulation of bone 
resorption is realised in the model via modulation of the 
rank— rankl— opg signalling pathway by the microscopic 
strain energy density ipbm- The production rate of rankl on 
pre-osteoblasts, [ s assumed to be enhanced during 

x 1 RANKL 3 a 



4' 



bin 



^bm(to) 



0, 



V'bm < ^bm(to) 
^bm > V'bmOo) 



(10) 



where k is a parameter quantifying the strength of the 
biomechanical transduction and ipbm^o) is tne steady-state 
value of the strain energy density. The steady state is 
assumed to be an initial homeostatic state of bone remod- 
elling, with no bone gain or loss. The increase in rankl 
production rate during mechanical disuse increases both 
and 71™ in Eqs (7) and (8), and consequently 

bcl,0C[j actjOCp 

leads to increased osteoclast generation (see Appendix A, 
Eqs. (18), (22)). 

Following Ref. [22], the biomechanical regulation of 
bone formation is realised in the model via modulation of 
pre-osteoblast proliferation by the microscopic strain en- 
ergy density ipbm- ^ n Eq- (5), pre-osteoblasts are generated 
both by differentiation from ob u s and by self-expansion 
through proliferation. The modulation of proliferation by 
the strain energy density is expressed by an 'activator' 
function U^n m , defined as: 



x act,OB p 
( 1 

2 



n 



act,OB p 



4' 



bm 



1 A 

2 + 2 VV'b m (t ) 

1 



4>bm < V'bmOo), 

1 ), ^bm(to) < ^bm < ^bm*, 
4>bm* ^ 4>bm, 

(ID 



where A is a constant parameter quantifying the strength 
of the biomechanical transduction and, i/>b m * = (1 + 
A _1 )i/j bm (t ) is the minimum value of the strain energy 
density for which n^™ B = 1. 



2.3.2 Geometrical feedback regulation 

The four regulatory functions g 0Bii (S v ), g 0Bp (S v ), g OCu (S v ) 
and g oc (S v ) in Eqs (5)-(8) include a geometrical feedback 
at various stages of osteoblast and osteoclast development. 
This enables us to distinguish two types of geometrical 
action. Indeed, the modulation of the bone cell develop- 
mental stages by the specific surface can be interpreted in 
the model in terms of modulation of the initiation of new 
remodelling events (bmu creation) and modulation of re- 
sorption and formation within existing bmus, as explained 
in the following: 

1. Initiation of new remodelling events. The initiation of 
a new remodelling event is a localised process that 
creates a new bmu. While the exact biochemical 
mechanisms that lead to the creation of a new bmu 
are poorly understood, it is believed that first steps in 
this process are the recruitment of pre-osteoclasts and 
pre-osteoblasts at the bone surface [3]. This recruit- 
ment is thought to be controlled by osteocytes sensing 
the local mechanical state of the bone matrix and 
communicating with progenitor cells in the marrow 
through the bone surface. The complex dependence 
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of bone surface availability on the recruitment of pre- 
osteoblasts and pre-osteoclasts is modelled by the 
geometrical regulation of cell differentiation exerted 
by the functions g OBu (S v ) and g oc (S v ). Therefore, 
the geometrical regulation by g 0B (s v ) and g 0Cu (S v ) 
models the influence of bone surface availability on 
the initiation of new remodelling events, and so on 
the number of bmus in the representative volume 
element, which in turn determines the bone turnover 
rate [2] . This type of geometrical regulation of bone 
remodelling is similar to the geometrical regulation 
of the 'activation frequency' of bmus used by Hazel- 
wood etal. [21]. 

2. Modulation of resorption and formation within exist- 
ing bmus. Active osteoclasts can only resorb bone 
from the bone surface. Similarly, active osteoblasts 
are only observed to deposit new bone at the bone 
surface. The maturation of pre-osteoblasts and pre- 
osteoclasts into active cells thus depends on bone sur- 
face availability. This dependence is modelled by the 
geometrical regulation of cell differentiation exerted 
by the functions g OB (S v ) and g oc (S v ). Therefore, 
these functions dictate how many active osteoclasts 
and active osteoblasts can form in bmus that are 
already remodelling bone (i.e., which contain already 
pre-osteoblasts and pre-osteoclasts), and so how much 
bone is resorbed and formed in existing bmus, which in 
turn determines bone balance. Note that while bone 
surface availability is necessary for cell activation, it 
is not sufficient. For instance, osteoclasts can only 
become active if, in addition, their receptor rank is 
activated by the ligand rankl. 

Due to the different ways of action of bone surface on cell 
differentiation explained above, the regulatory functions 
g B u (Sv), g OBp (Sy), g 0Cu (S v ) and g 0Cp (Sv) may take differ- 
ent forms and can be complicated functions of the specific 
surface S v . We take here a phenomenological approach and 
assume that each of these functions can be represented by 
a power law of S v : 

f S v \k t 

g,( S v) = I c 77 S j ' With k i - °> i = 0B u> 0B P> 0C u> 0C P- 

(12) 

In Eq. (12), S v [t ) denotes the specific surface in the 
bone remodelling steady state, which is assumed to be 
homeostatic (no bone gain or loss, but a steady bone 
turnover). The normalisation of S v by its steady-state value 
S v (t ) ensures that in the steady state, g, = 1 Vi, and 
so that the steady state of the model is consistent with 
previous models of bone remodelling without geometrical 
regulation [22,23]. The benefit of using a power-law 
function of S v in Eq. (12) is that geometrical feedback 
can be "switched off" by choosing fc ; = Vi in Eq. (12). 
In this situation, not only the steady state, but also the 
dynamical behaviour of the model of Ref. [22] (including 
mechanical feedback) is retrieved. Mechanical feedback 
is "switched off" by choosing k: = 0, A = in Eqs (10) 
and (11). Note that when both geometrical and mechanical 
feedbacks are "switched off", the bone cell population 



model of Refs [23, 24] is retrieved, except for a remaining 
pre-osteoblast proliferation term in Eq. (5) that was not 
accounted for previously. 5 

To reveal in which ways the morphology of the mi- 
crostructure of bone may influence bone remodelling, we 
will investigate the effects of several combinations of g 0B 
goB > Soc an d goc i n Section 3 and determine combinations 
that lead to physiologically meaningful results. 

2.4 Changes in porosity and bone matrix 
fraction due to cell activity 

The activity of osteoclasts and osteoblasts leads to the 
removal and deposition of new bone. This activity modifies 
the volume fraction of bone matrix in the tissue. Os- 
teoblasts deposit osteoid, a collagen-rich substance which 
later mineralises into new bone. Primary mineralisation 
of osteoid is relatively fast: 70% of the maximum mineral 
density is reached within a few days in humans [2]. Given 
the much larger time spans involved in the remodelling 
processes, it is fair to model the osteoblasts "instanta- 
neously" depositing "fully" mineralised new bone matrix, 
as was assumed in Refs [22-24]. We further assume that 
the resorption rate of bone matrix fc res by an individual 
active osteoclast (in volume per unit time) is constant, and 
that the rate of new bone matrix deposition fc form by an 
individual osteoblast (in volume per unit time) is constant. 
The evolution of the vascular porosity and bone matrix 
volume fraction are thus given by 

d _ d 

T^/vas — — j^ bm ~~ — ' c form OB a + Kes oc a- (13) 

2.5 Comparison with the model of geomet- 
rical regulation of bone remodelling by 
Martin 

It is informative at this point to compare our formulation 
of geometrical regulation of bone remodelling with that 
proposed by Martin in Ref. [13]. The evolution of the 
vascular porosity proposed by Martin is, in our notations 
(see [13, Eq. (67)]): 6 

d 

Tlfvas, Martin — — C^form ^ob. ^-ob. — ^-res ^oc, ^oc.) (14) 

dt a a a a 

where 5 OBa and 5 0Ca are the active osteoblast and active os- 
teoclast surface densities at an active bone site (number of 
cells per unit surface), and A OBa and A OCa are the fractions of 
the total available bone surface in which there is osteoblas- 
tic and osteoclastic activity. The quantity fc form 5 0Ba A 0Ba S v 
is the formation rate of bone matrix and corresponds to 
fc form OB a in our model. The quantity fc res 5 OCa A OCa S v is the 
resorption rate of bone matrix and corresponds to A: res oc a 

5 Compared to Refs [23,24], the differentiation rate of ob u s to OB p s 
is reduced accordingly to ensure that the model converges to the same 
steady state. 

6 In Ref. [13], Martin uses the symbol S% to denote only the fraction of 
3D specific surface that is capable of remodelling. This excludes surfaces 
of the lacunae-canaliculi system [17]. In this paper, we do not consider 
surfaces of the lacunae-canaliculi system, and so in Ref. [13] is equal 
to our definition of S v . 
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in our model. To model an osteoporotic pathological con- 
dition, Martin assumes a constant bone imbalance between 
resorption and formation in remodelling bmus. This imbal- 
ance is modelled by setting (fc form 5 0Ba A 0Ba -fc res 5 0Ca A OCa ) = 
—2 p,m/year in Eq. (14), meaning that in average, a 2 Lim- 
thick layer of the bone surface is resorbed each year. The 
results presented in Figure 1 were obtained by Martin in 
Ref. [13] by integrating Eq. (14) with this constant imbal- 
ance and with either S v = 1 (no geometrical feedback) or 
using the function S v (/ vas ) presented in Eq. (4) to include 
the effect of geometrical feedback. 

A limitation of the formulation by Martin is that the de- 
pendence of the rate of change of porosity upon S v does not 
account for microscopic effects of bone surface availability 
on the recruitment and development of bone cells. Indeed, 
in Martin's formulation, the complexity of osteoblast and 
osteoclast recruitment and development at a remodelling 
bone surface is embodied by the fractions A 0Ba and A 0Ca , 
but these fractions are assumed to be independent of S v . 
Martin and others have interpreted Eq. (14) as the influ- 
ence of geometrical regulation on the activation frequency 
of bmus [2,13]. Indeed, the densities of active osteoblasts 
and osteoclasts in Martin's model are proportional to the 
specific surface. Comparing Eqs (13) and (14), one sees 
that OB a = A OBa S v 5 OBa and oc a = A 0Ca S v 5 0cv Thus, an 
increased specific surface represents increased numbers of 
remodelling OB a s and oc a s, and so an increased number of 

BMUS. 

The strength of our model is to consider several stages in 
the recruitment and development of osteoblasts and osteo- 
clasts explicitly. This enables us to include a more detailed 
geometrical feedback acting directly on these different 
stages of the remodelling sequence. In our model, even 
steady-state values of OB a s and oc a s are complex functions 
of the specific surface. This complexity represent the im- 
plicit dependence of A 0Ba and A OCa upon S v not accounted 
for by Martin. Additionally, no mechanical feedback has 
been considered by Martin. In Ref. [22], Schemer etal. 
have discussed the problem of neglecting mechanical feed- 
back in models of bone remodelling applied to catabolic 
bone diseases (such as osteoporosis), which leads to a 
unlimited bone loss which is physiologically not observed. 
This behaviour can be seen in Figure 1 for cortical bone 
where the vascular porosity continues to increase (both 
with and without geometrical regulation). It is known 
from bone mineral density measurements in osteoporotic 
patients that the increase in bone porosity in osteoporosis 
slows down with time. The porosity eventually reaches an 
upper bound whose value depends on the patient. We note 
here that the long-term value of the vascular porosity / vas 
reached with mechanical feedback in our model depends in 
particular on the initial value of / vas (see Section 3). 

3 Numerical Results and discussion 

Having incorporated biochemical, biomechanical and ge- 
ometrical regulation of bone remodelling we are now in 
a position to investigate the effects of various regulatory 
parameters on changes in bone cell numbers, vascular 
porosity and bone stiffness. To compare our model with 



the model suggested by Martin in Ref. [13], we simulate an 
underlying osteoporotic condition as a perturbation from 
the original (homeostatic) steady state situation. 

Osteoporosis is a bone disease that leads to an increase 
in porosity in both cortical and trabecular bone [18]. This 
increase in porosity generates a progressive reduction in 
bone stiffness and a higher fracture risk. To simulate an 
osteoporotic condition in our model, we perturb the home- 
ostatic steady state (a state with no bone loss or gain and a 
constant bone turnover rate) by increasing the rankl/opg 
ratio. This can be achieved in our model by prescribing 
an excess of pth concentration (see Refs. [23, 25]). 7 The 
osteoporotic condition is assumed to develop instantly from 
an initial homeostatic state at time t (corresponding to 
a middle biological age). To obtain a steady increase in 
/ vas of 0.01/year without geometrical regulation, as has 
been assumed by Martin (see Figure 1), a continuous pth 
administration rate P PTH (t) = 500 pM/day has been applied 
at all times t > t . We denote by At 0P = t — t the time 
elapsed since the onset of the osteoporotic condition. The 
evolution of the system is followed for a period of time of 
20 years from the onset of osteoporosis, i.e., < At 0P < 
20 years. The initial values for vascular porosity of cortical 
and trabecular bone have been chosen as / v ™ rt ' = 0.05 and 
fvas h =0.75. For the simulations including biomechanical 
regulation, a uniaxial compressive stress S = S um e 3 ® e 3 
with E um = — 30 MPa has been assumed to exert on the 
representative volume element. 

3.1 Simulation of osteoporosis: Evolution of 
bone porosity and bone stiffness proper- 
ties 

Figure 5 shows the simulated evolution of bone cell den- 
sities, vascular porosity (/ vas ) and selected components of 
the bone stiffness matrix (i.e., axial stiffness E 3 and shear 
modulus G 12 ) in osteoporosis considering no geometrical 
feedback and no mechanical feedback. The evolution of 
bone cell densities exhibits a short transient for a period 
of as 30 days after the onset of osteoporosis (Figure 5a). 
After this initial transient, oc a s and OB a s reach a new steady 
state, in which resorption and formation are imbalanced 
and lead to osteoporotic bone loss. Note that a short 
time interval has been chosen in Figure 5a to demonstrate 
the short transient cellular response, while a larger time 
interval of 20 years is chosen to follow the evolution of 
the vascular porosity, axial stiffness and shear modulus. 
The progressive increase in vascular porosity is shown 

7 The physiological effect of pth administration for bone remodelling 
is complex and in particular, depends on the time course of the admin- 
istration [34]. Continuous pth administration (infusion) leads to bone 
loss with increased turnover. However, intermittent pth administration 
(daily injections) leads to bone gain [34, 35] . Only the 'continuous' action 
of pth is represented in our model [23, 25]. We have shown previously 
that increasing pth leads to an increase in the rankl/opg ratio, and so to 
bone loss with higher turnover rate compared to the homeostatic bone 
remodelling state [24,28]. Increasing pth concentration thus consists 
in an adequate representation of osteoporosis in the present model. We 
also note that the effect of pth in bone remodelling is known to interact 
synergistically with mechanical loading [36, 37]. Consistently with the 
interpretation of increasing pth as a representation of osteoporosis in our 
model, we do not consider this synergistic interaction here. 
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Figure 5 - Time evolutions for a simulated osteoporotic condition of (a) 
bone cell densities (normalised to the initial homeostatic steady state) and 
(b) the vascular porosities / vas , axial stiffnesses E 3 and shear modulus G 12 
in both a cortical bone tissue sample and a trabecular bone tissue sample 
under uniaxial compression. In this figure, no geometrical nor mechanical 
feedback is considered. 



in Figure 5b together with the associated reduction in 
bone stiffness. After 20 years of osteoporosis, the normal 
stiffness in the longitudinal direction of cortical bone, E 3 , 
has dropped by about 20%, whereas the shear modulus of 
cortical bone, G 12 , has dropped by about 40%. From these 
results it is clear that the evolution of vascular porosity 
drives the changes in bone stiffness differently for the 
different components of the stiffness tensor [28]. Generally 
speaking, an increase in / vas is always associated with a 
reduction of bone stiffness properties. For conciseness, in 
the following we will only present the effects of geometrical 
and biomechanical regulation on the vascular porosity / vas 
alone. 



3.2 Geometrical regulation — effect on indi- 
vidual stages of osteoblast and osteoclast 
developments 

We first present how geometrical regulation influences the 
osteoporotic increase in bone porosity / vas without account- 



ing for mechanical feedback, i.e., by setting k = pM/day 
and A = in Eqs (10), (11). Based on the four regulatory 
functions g 0Bu , g 0Bp , g OCy , and g 0Cp accounting for geomet- 
rical regulation, we first performed 2 4 = 16 simulations 
'switching on' or 'off each of these regulatory functions to 
investigate all different possible combinations. However, it 
turns out that among all possible simulations there are only 
three patterns which reflect the bone systems behaviour. 
These patterns can best be studied by looking first at the 
cases in which only a single regulatory function is 'switched 
on' while all others are 'switched off (i.e., identically 
set to one) (Figure 6). We will present a selection of 
combinations of geometrical regulation on several stages of 
osteoblast and osteoclast development in Figures 7 and 8. 

Figure 6a shows the influence of g 0Bu (S v ) on the increase 
of vascular porosity in osteoporosis compared to the case 
of no geometrical and no biomechanical regulations in 
both cortical and trabecular bone. The time evolution of 
the porosity / vas (t) in Figure 6a clearly shows that while 
g OBu [S v ) reduces bone loss in cortical bone it accelerates 
bone loss in trabecular bone. The strength of this reg- 
ulation depends on the exponent k 0Bp of the normalised 
geometrical regulatory function in Eq. (12). Exponents 
in the range 0.3 < fc 0B < 1 exhibit a relatively strong 
regulatory effect on the evolution of the disease, while 
exponents in the range < fc 0B < 0.2 only exhibit a 
moderate effect. Interestingly, the mechanism of action 
of g 0B is opposite to the geometrical regulation obtained 
by Martin [13] (see Figure 1). In cortical bone, the 
osteoporotic increase in / vas induces an increase in S v (see 
Figure 4), and so an increase of g 0Bu (S v ). This increases in 
turn the generation of osteoblasts, which has an anabolic 
(i.e., bone forming) effect and stabilises the osteoporotic 
loss of bone [23, 25]. By contrast, in trabecular bone the 
osteoporotic increase in / vas induces a decrease in S v (see 
Figure 4), and so a decrease of g 0B (S v ). This decreases 
in turn the generation of osteoblasts, which accelerates the 
osteoporotic loss of bone. 

While this scenario can be understood from a theoretical 
point of view, it has to be emphasised that it is probably 
physiologically unrealistic. As argued in Section 2.3.2, the 
geometrical regulation of ob u differentiation is associated 
with the creation of a new remodelling event, i.e., of 
a new bmu. However, both ob u differentiation and oc u 
differentiation are activated in such an event. Here, the 
consideration of a geometric regulation on ob u differenti- 
ation alone represents the initiation of a formation event 
and is not associated with a joint initiation of a resorption 
event. To represent the geometrical regulation of bmu 
creation, both g 0Bu (S v ) and g 0Cu (S v ) should act together 
(see Figures 7 and 8). 

Figure 6b shows the influence of g 0B (S v ) on the increase 
of / vas in osteoporosis. In this figure, / vas is only weakly 
affected by the geometrical feedback in both cortical and 
trabecular bone. The effect of bone surface availability on 
the differentiation of pre-osteoblasts into active osteoblasts 
modelled by g OB (S v ) does not seem to affect the bone 
remodelling balance significantly in the first 10 years. In 
fact, neither the population of active osteoblasts nor the 
population of active osteoclasts is affected significantly, 
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Figure 6 - Influence of geometrical regulation on the increase of vascular porosity in osteoporosis. In this figure, geometrical regulation is included at a 
single stage of osteoblast or osteoclast development only. The numerical simulations are initiated from either cortical bone (lower curves) or trabecular 
bone (upper curves). The effect of several strengths k of the geometrical regulations are shown (corresponding to the exponent of the power law in 
Eq. (12)). The case k = corresponds to no geometrical feedback, (a) Geometrical regulation of ob u differentiation only; (b) geometrical regulation of 
OB p differentiation only; (c) geometrical regulation of oc u differentiation only; (d) geometrical regulation of oc p differentiation only. 



while the population of pre-osteoblasts is decreased in cor- 
tical bone and increased in trabecular bone (not shown). 
The increased differentiation rate of pre-osteoblasts into 
active osteoblasts in cortical bone (due to the increase 
in the available bone surface) depletes the population of 
pre-osteoblasts. The population of active osteoblasts is 
thus derived from a smaller pool of pre-osteoblasts that 
are differentiating faster, and therefore stays relatively 
constant. The opposite effect occurs in trabecular bone, 
i.e., a larger pool of pre-osteoblasts is created but they are 
differentiating into active osteoblasts more slowly. 

Figure 6c shows the influence of g QCu (S v ) on the increase 
of / vas in osteoporosis. The behaviour of / vas obtained for 
this type of geometrical regulation is in general agreement 
with the geometrical regulation obtained by Martin [13] 
(see Figure 1), i.e., geometrical regulation leads to in- 
creased bone loss in cortical bone and decreased bone 
loss in trabecular bone. In cortical bone, the increase 
in g oc (S v ) due to the osteoporotic condition leads to 



an increase in osteoclastogenesis, which by biochemical 
coupling also increases (although to a lesser extent) the 
population of osteoblasts. This results in a high turnover 
rate with a catabolic bias, i.e., a high rate of bone loss. 
In trabecular bone, the situation is reversed since g OCu (S v ) 
decreases with the evolution of the osteoporotic condition. 
A state with lower turnover rate and a reduction in resorp- 
tion is reached. The low turnover rate explains the reduced 
amplitude of the response in trabecular bone compared to 
that seen in cortical bone. 

The strength of the geometrical regulation is determined 
by the value of fc OCu and is seen to strongly affect bone 
balance. In particular for the case of strong modulation, 
i.e., fc oc = 1, a rapid acceleration of cortical bone loss 
occurs over the first 10 years of osteoporosis. At nearly 
10 years, a transition is taking place leading to a reduction 
in the rate of bone loss. This behaviour is due to the 
fact that at this time, the bone porosity reaches the critical 
value / vas * ^ 0.37 at which the specific surface is maximum 



10 



(see Figure 3). In the first 10 years of osteoporosis, / vas 
is in the ascending branch of S v (/ vas ) and so g OCu (S v ) in- 
creases, while after reaching the maximum specific surface, 
goc u ($v) decreases, which leads to a reduction in the rate 
of bone loss. Over 20 years, the original cortical bone 
has been resorbed enough to reach trabecular porosities. 
Such a strong effect of bone remodelling on bone porosity 
is normally not physiologically seen. However, it is possible 
that in osteoporosis, a locally strong geometrical feedback 
may help initiate or accentuate the observed 'trabeculari- 
sation' of bone, by which cortical bone is progressively lost 
and transformed into trabecular bone at the endocortical 
wall, leading to cortical wall thinning and expansion of 
the medullary cavity [19,38-40]. Indeed, in bone, the 
endocortical wall exhibits the highest specific surface and 
is known to be highly remodelling. It can be expected that 
geometrical regulation plays a particularly significant role 
in this region of bone. 

Figure 6d shows the influence of g oc (S v ) on the increase 
of / vas in osteoporosis. This influence is qualitatively similar 
to that of g oc shown in Figure 6a, i.e., a reduction in 
bone loss in cortical bone and an increase in bone loss 
in trabecular bone. However, the geometrical regulation 
modelled by g oc is less pronounced than the geometrical 
regulation modelled by g 0Bu . 

3.3 Geometrical regulation — combined ef- 
fect on several stages of osteoblast and 
osteoclast developments 

As mentioned previously, a geometrical regulation of the 
creation of new remodelling events (i.e. new bmus) should 
involve a regulation of both the recruitment of osteoclasts 
and the recruitment of osteoblasts (see Section 2.3.2). 
In Figure 7a, we show the combined influence of both 
g 0Cu (S v ) and g 0B (S v ) on the increase of / vas in osteoporo- 
sis. The exponents fc oc and k 0% measuring the strength 
of the geometrical regulation are assumed identical. By 
comparing with the individual influences of g OBu (S v ) and 
g 0Cu (S v ) in Figure 7a and 7c, one sees that the geomet- 
rical regulation of osteoblast recruitment overrides that 
of osteoclast recruitment. As a consequence, the overall 
behaviour is opposite to that obtained by Martin [13] for 
the geometrical regulation of the creation of new bmus (see 
Figure 1). 

By contrast, in Figure 7b, the geometrical regulation 
of the last stage of osteoclast differentiation by g oc (S v ) 
and the geometrical regulation by g 0B (S v ) of the last 
stage of osteoblast differentiation (modelling a regulation 
of the activation of cells within existing bmus) seem to 
compensate each other and to result in an evolution of 
the porosity that is almost unaffected by geometrical feed- 
back. The geometrical regulations assumed in Figures 7a 
and 7b represent different natures of the influence of bone 
surface availability on bone remodelling (see Section 2). 
Therefore, our simulations suggest that for the simulated 
evolution of bone porosity in osteoporosis, the influence of 
surface availability is significantly stronger on the creation 
of new remodelling events (new bmus) (Figure 7a) than on 
the activation of bone cells within already active remod- 
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Figure 7 - Combined influence of several geometrical regulations on the 
increase of vascular porosity in osteoporosis. The numerical simulations 
are initiated from either cortical bone (lower curves) or trabecular bone 
(upper curves). The effect of several strengths k of the geometrical 
regulations are shown (corresponding to the exponent of the power law 
in Eq. (12)). The case k = corresponds to no geometrical feedback, 
(a) Joint geometrical regulation of both ob u differentiation and oc u 
differentiation; (b) joint geometrical regulation of both OB p differentiation 
and oc p differentiation. 

elling sites (within existing bmus) . 

In Figure 8, we show that a similar influence of geo- 
metrical regulation of bmu creation as that obtained by 
Martin [13] can be retrieved in our model by modifying 
the relative strengths of the regulatory functions g 0Bu (Sv) 
and g 0Cu [S v ) via the exponents fc OBi and fc 0Cu . The vas- 
cular porosity can exhibit a wide range of behaviours, 
interpolating between the situation of Figure 6a (fc 0Bi = 
1, fc 0Cu = 0) and the situation of Figure 6c (fc 0Bu = 0, 
fc 0Cu = 1). We conclude that geometrical feedback has the 
potential to significantly influence the evolution of bone 
diseases. However, complex biochemical coupling between 
osteoblasts and osteoclasts makes it difficult to predict the 
relative strength of the influence of geometrical regulation 
on osteoblast and osteoclast developments. 

Finally, we note that because the curve S v (/ vas ) exhibits a 
maximum at/ vas * 0.37 (see Figure 3), geometrical feed- 
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Figure 8 - Influence of the geometrical regulation of bmu creation on the 
increase of vascular porosity in osteoporosis. Different strengths for the 
geometrical regulation of ob u differentiation and oc u differentiation are 
investigated via the exponents k OBa and fc OCu of the regulatory functions 
Sob u (SvO and g 0Cu (S v ): k 0H = 0,k OCu = (no geometrical regulation, 
thin dashed lines); fc 0B = l,fc oc = 0.1 (thick dashed lines) and ic 0B = 
0.1,fc OCu = 1 (thick solid line). 

back always induces opposite behaviours for the evolution 
of vascular porosity in cortical bone (for which / vas < / vas *) 
and in trabecular bone (for which / vas > / vas *) . 

3.4 Coupled geometrical and mechanical 
regulations 

Figure 9 shows the effect of adding a mechanical regulation 
of bone remodelling for the evolution of the vascular poros- 
ity. The geometrical regulation considered in Figure 9 (see 
dotted line) is assumed to represent the influence of bone 
surface availability for bmu creation as in Figure 8, solid 
line. Two strenghts of the biomechanical transduction are 
illustrated (i.e., X = 0.1 and X = 0.5 in Eq. (11)). Figure 9 
suggests that mechanical feedback has the potential to pro- 
gressively override both the evolution of osteoporosis and 
the influence of geometrical regulation. Indeed, mechani- 
cal feedback is seen to stabilise bone balance, irrespective 
of whether geometrical feedback stabilises or destabilises 
bone balance in osteoporosis. Such a stabilisation of bone 
loss is clinically observed in osteoporotic patients. 

Similarly to the simulations by Scheiner etal. [22], 
mechanical feedback counteracts bone loss in osteoporosis 
both in cortical and trabecular bone, on a time scale that 
depends on the strength of mechanical regulation (i.e., on 
the parameter A). Interestingly, in our simulations with 
geometrical feedback, the strength X of the mechanical 
regulation also has an influence on the steady-state value 
of the porosity attained. This is to be contrasted to the 
situation without geometrical feedback in which no such 
influence was observed [22]. In fact, without geomet- 
rical feedback, the vascular porosity f vas enters the right 
hand side of the governing equations (5)-(8), (13) only 
implicitly via the strain energy density xp bm , see Eq. (9). 
Consequently, the steady-state value of _f vas is uniquely 
determined by the steady-state value of the strain energy 
density, ipt, m [t Q ), and the macroscopic loading S, which are 
themselves independent of the strength X of biomechanical 
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Figure 9 - Influence of coupled geometrical and mechanical regulation of 
bone remodelling on the increase of vascular porosity in osteoporosis. The 
evolution of / vas without geometrical regulation and without mechanical 
regulation is shown as dashed lines. The inclusion of geometrical 
regulation without mechanical regulation is shown as dotted lines. This 
geometrical regulation is assumed identical to the solid line in Figure 8, 
i.e., fc OBu = 0.1, fc 0Cu = 1, k 0Bp = fc OCp = 0. Two cases of coupled 
geometrical and mechanical regulations are shown as solid lines for two 
strengths A of the biomechanical regulation: A = 0.1 (thin solid lines) and 
X = 0.5 (thick solid lines). 

regulation [22]. With geometrical feedback, additional 
dependences upon the vascular porosity / vas are added in 
the governing equations for the bone cells, Eqs (5)-(8). 
Consequently, the steady-state value of / vas now depends 
in addition on the biochemical and cellular state of the 
system, which depends in turn on the strength X of the 
biomechanical regulation of the bone cells. 

From the above discussion, it is clear that the coupling 
between the geometrical and mechanical feedbacks is effec- 
tively mediated by the biochemical and cellular state of the 
system. In osteoporotic patients, both a geometrical and a 
mechanical feedback are present at the same time, as well 
as biochemical and hormonal dysregulations underlying 
the establishment of osteoporosis. The interdependence 
between the biochemical and hormonal dysregulations in 
osteoporosis, geometrical feedback and mechanical feed- 
back, is not trivial to elucidate. Physiologically, it is ex- 
pected that all these influences play a role for the evolution 
of bone vascular porosity. Our mathematical model is a first 
attempt to integrate these influences and help understand 
their contribution for clinically observed changes of bone 
porosities in osteoporotic patients. 

4 Conclusions 

In this paper we developed a computational model of bone 
remodelling that takes into account biochemical, biome- 
chanical and geometrical regulations of bone cells. The 
biochemical regulation of the bone cells is based on the 
model developed in Refs [23, 24]. The biomechanical 
regulation of the bone cells is based on the model devel- 
oped in Ref. [22], in which a continuum micromechanical 
approach is used to consistently link bone cell responses 
with mechanical properties of bone. The new contribution 
of this paper is the inclusion of a geometrical regulation 
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of bone cells in the model. Geometrical feedbacks were 
included at several developmental stages of osteoblasts 
and osteoclasts to represent the influence of microscopic 
bone surface availability for various bone microstructures. 
We investigated the influence of a geometrical regulation 
of bone cells in bone remodelling both without and with 
consideration of biomechanical regulation for a simulated 
osteoporotic condition. From the numerical simulations, 
we identified the following actions of geometrical and 
mechanical regulation on bone remodelling: 
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• Geometrical regulation of bone remodelling may play 
an important role for the initiation of new bmus as 
described by the combined effect of the geometri- 
cal regulatory functions g 0Bu and g oc acting on ob u 
differentiation and oc u differentiation; in particular 
geometrical regulation of oc u s via g 0Cu seems to be 
most important to retrieve similar evolutions of bone 
porosities in osteoporosis as obtained by Martin [13]; 

• Geometrical regulation of bmu creation affects cortical 
and trabecular bone in opposite ways in osteoporosis: 
while bone resorption is enhanced in cortical bone 
due to the increase in specific surface with increasing 
porosity, bone resorption is slowed down in trabecular 
bone due to the decrease in specific surface with 
increasing porosity; 

• Geometrical regulation of the activation of osteoblasts 
and osteoclasts in existing bmus seems to play a sec- 
ondary role for the evolution of osteoporosis. While 
the specific surface can influence the differentiation 
of bone precursor cells into active resorbing/forming 
cells, no significant influence on bone porosity was 
observed in our simulations. 

• Our simulations suggest that geometrical regulation 
may play a role in osteoporosis for the initiation and/ 
or accentuation of the observed 'trabecularisation' of 
bone at the endocortical wall. At the endocortical wall, 
the specific surface and bone turnover are high and so 
the effects of a geometrical feedback can be expected 
to be significant; 

• Our simulations of coupled geometrical and mechan- 
ical regulations suggest that the stabilisation of bone 
loss observed clinically in osteoporotic patients is 
probably accelerated by geometrical feedback in tra- 
becular bone, but is probably slowed down by geo- 
metrical feedback in cortical bone. 



Appendix A Model description 

In this appendix, we complete the mathematical descrip- 
tion of the model and give the values of the parameters. 
The nomenclature used in the paper is split into three 
tables. Table 1 lists dynamics quantities involved in the 
governing equations of the bone cells and bone porosity, 
Eqs (5)-(8), (13). Table 2 lists all the parameters relevant 
to the biochemical regulation of the model. Table 3 lists 
quantities involved in the biomechanical regulation of the 
model. 



Concentrations of the biochemical signalling molecules 

As in Ref. [22-24,26], the concentration of the biochemical 
signalling molecules are governed by rate equations based 
on mass action kinetics. Ligand-receptor binding reactions 
occur on a time scale much faster than the characteristic 
times of cellular response (such as differentiation, apop- 
tosis). The rate equations for the biochemical signalling 
molecules can therefore be taken in their steady state (see 
Refs. [23,26] for details). This gives: 



TGF P (t) = nt°" e fc res oc a (t)/D. 



bone 7 
'tgf[3 ^res*- 

RANK(t) = iV o 7 I< OC p (t), 



TGF[3 



OPG(t) : 



£o°B P a G 0B aUK™, 0B 



^OB a (07r^V OPG max + Do 



RANKL e ff 



6 +p^b 

r'RANKL ' r ran 



RANKL 4" "^RANKL RANKL eff 



RANKL(t) = 

1 4" "^[rANKL-OPg] OPG + -K[ RANK L-RANK 

PTH(t)= [P PTH (t) + /3 PTH ]/D 

PTH ' 



RANK 



(15) 
(16) 

(17) 



(18) 



(19) 



In Eq. (18), RANKL eff is the maximum concentration of 
rankl (also referred to as effective carrying capacity), 
which is regulated by the parathyroid hormone pth: 



_ . - TT , T atRANKL „_ —PTH 

RANKL eff = K« 0B P ^act.o 



(20) 



Both mechanical and geometrical feedbacks are im- 
portant to account for in our model of bone remod- 
elling. Mechanical feedback enables the local porosity 
of bone tissue / vas to stabilise to a well-defined value 
within [0,1]. Geometrical feedback enables this value 
to be not only determined by the external loading of 
the tissue, but also by the biochemical and cellular 
state of the system, as would be expected physiolog- 
ically 



where jv RANKL is the maximum number of membrane-bound 

OB p 

rankl molecules that can be expressed on a single pre- 
osteoblast (see Refs. [23, 24] for more details). In this 
work, the concentration of mcsf, and so the quantity 



r MCSF 

act,oc u 3 



, are assumed constant (see below and Table 2). 



The additional production rate of pth in Eq. (19), P PTH (t), 
is used to increase the normal systemic levels of pth 
to simulate an osteoporotic condition (see comments in 
Section 3). 
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Activator and repressor functions The regulation of the 
bone cell behaviours (such as differentiation, apoptosis, 
expression rate of a ligand) by the biochemical signalling 
molecules is modelled in Eqs (5)-(8) by so-called "ac- 
tivator" functions and "repressor" functions rcj ' 
where L denotes the signalling molecule (ligand) and A 
the signalled cell. These activator and repressor functions 
represent the strength of the response of the cell to the 
signal mediated by the ligand and are assumed to be given 
by: 



n 



act,A 



rep,A " 



act,A 



(21) 



where is the dissociation binding constant between the 
ligand and its receptor on the cell. The quantity 7Tg CtA 
represents the fraction of the receptors on the cell A that are 
bound to a ligand L (see Refs. [23,25] for more details). 



For example, the activator functions n 
defined as: 



RANKL 

act,oc.. 



and n 



act,oc p 



are 



„ RANKL „ RANKL 

71 =71 

act,oc u act,oc p 



RANKL 



RANKL + 



where k*f KL is the dissociation binding constant between 
rankl and the rank receptor on oc u s and oc p s, and rankl is 
the free (unbound) rankl concentration given by Eq. (18). 

Steady-state values The geometrical and biomechanical 
regulations of the bone cells is normalised by the steady- 
state values. These values depend in particular on the 
initial bone porosity, and so are calculated prior to evolving 
the system. 

Table 1 - Dynamic quantities in the governing equations, Eqs (5)-(8), 
(13) 

SymbolUnit Description 

oc p pM density of pre-osteoclasts 

oc a pM density of active osteoclats 

ob p pM density of pre-osteoblasts 

OB a pM density of active osteoblasts 

tgf(3 pM concentration of transforming growth factor p 
rank pM concentration of receptor-activator nuclear factor kB 
rankl pM concentration of receptor-activator nuclear factor jcB 
ligand 

opg pM concentration of osteoprotegerin 
pth pM concentration of parathyroid hormone 



fvas 



volume fraction of vascular pores 
specific surface 



4\ 



MPa microscopic strain energy density of the bone matrix 



Symbol 



D 
D 

D n 



(22) n 



Table 2 - Biochemical parameters 



Value Description 



oc u 1 • 10" 3 pM 

ob u 1 • 10" 3 pM 

fc res 200 pM^day -1 

k form 40 pM-May" 1 



density of uncommitted osteoclast 
progenitors 

density of uncommitted osteoblast 
progenitors 

daily volume of bone matrix re- 
sorbed per osteoclast 
daily volume of bone matrix formed 
per osteoblast 



^. RANKL 
OC 

7.TGFp 

7,TGFp 
K OB„ 



1LPTH 

OB,act 



J^PTH 

OB,rep 



^[rankl-rank] 



^[rankl-opg] 



oc u — > oc p differentiation rate pa- 
rameter 

oc p — > oc a differentiation rate pa- 
rameter 

oc a apoptosis rate parameter 
ob u — > OB p differentiation rate pa- 
rameter 

0.166/day ob p — > OB a differentiation rate pa- 
rameter 

0.021/day OB p proliferation rate parameter 
0. Ill/day OB a apoptosis rate 



4.2/day 

2.1/day 

5.65/day 
0.7/day 



0.5 
5.68 pM 
5.63- 10" 4 pM 
5.63- 10" 4 pM 
1.75- 10" 4 pM 
150 pM 
0.222 pM 
0.034/pM 
0.001/pM 



for 



constant for 



constant for 



(>/■ 

ft: 



OOPG 

' °»a 

RANKL 
OB p 
jyRANK 
OC p 

OPG mav 



n. 



bone 

TGFp 



TGFp 



1.68 • 10 2 pM/day 
250pM/day 
500pM/day 

1.62- 10 8 /day 
2.7- 10 6 
1-10 4 
2- 10 8 pM 

0.01 pM 

2/day 
10/day 
0.35/day 
86/day 



value of the activator function of 
mcsf for oc u — > OCp differentiation 
dissociation binding constant for 
rankl binding on oc u and oc p 
dissociation binding constant 
tgf[3 binding on oc a 
dissociation binding 
tgfP binding on ob u 
dissociation binding 
tgfP binding on OB p 
dissociation binding constant for pth 
binding on ob (in tt™ ob ) 
dissociation binding constant for pth 
binding on ob (in tt™ ob ) 
association binding constant for 
rankl and rank 

association binding constant for 
rankl and OPG 
production rate of rankl 
production rate of systemic pth 
continous administration rate of pth 
to model osteoporosis 
production rate of opg per ob 3 
maximum number of rankl per ob p 
number of rank receptors per oc p 
opg density at which endogeneous 
production stops 

density of tgfP stored in the bone 
matrix 

degradation rate of tgf|3 
degradation rate of rankl 
degradation rate of opg 
degradation rate of pth 
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Table 3 - Biomechanical quantities 



Symbol Value Description 

S S um e 3 <g> e 3 , macroscopic stress tensor 

S uni _ _ 30 MPa 

e bm microscopic strain tensor of the bone 

matrix 

%p bm microscopic strain energy density of the 

bone matrix 

k 1 • 10 5 pM/day strength of biomechanical transduction 

of bone resorption 
X 0.1 or 0.5 strength of biomechnical transduction of 

bone formation 
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